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Atomic clouds prepared in 'timed Dicke' states, i.e. states where the phase of the oscillating atomic 
dipole moments linearly varies along one direction of space, are efficient sources of superradiant light 
emission [M.O. Scully et al., Phys. Rev. Lett. 96, 010501 (2006)]. Here, we show that, in contrast 
to previous assertions, timed Dicke states are not the states automatically generated by incident 
laser light. In reality, the atoms act back on the driving field because of the finite refraction of the 
cloud. This leads to non-uniform phase shifts which, at higher optical densities, dramatically alter 
the cooperative scattering properties, as we show by explicit calculation of macroscopic observables, 
such as the radiation pressure force. 
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I. INTRODUCTION 



The fact that ensembles of point-like scatterers respond collectively to incident radiation is well-known since the 
i— i' seminal paper of Dicke [lj. The collective phenomenon, termed superradiance, has been the topic of a huge amount 
\ of theoretical and experimental work. However, the question about the exact nature of the state generated by the 
. ' radiation travelling through the ensemble has only been raised very recently. Scully and coworkers pointed out Q , 
that the dipole moments of scatterers distributed along the incident beam's optical axis are excited in phase with 
beam's propagating phase front. The resulting 'timed Dicke' state emits light predominantly into forward direction, 
^ ' provided the ensemble size is large compared to the radiation wavelength [3j . 

This simple picture, however, only holds at low optical densities. At high optical densities and spherically symmetric 
ensembles, the refraction index of the scattering medium delays the propagation of the pump beam and distorts its 
phase front. This distortion can have a huge impact on any macroscopic observable of the system, such as the 
angular distribution of the scattered radiation, the phase-front of the transmitted beam or the force acting on the 
Qh' center-of-mass of the ensemble. 

The aim of this paper is to calculate the correct state generated in an ensemble of two-level systems, e.g. an 
' atomic cloud, by interaction with a weak laser beam and how this state cooperatively scatters the incident light. 
, Cooperative scattering by many atoms has been studied extensively in the past, both classically and quantum- 
t — 1 mechanically. Classically, scattering at extended objects is described as Mie scattering, showing resonances induced 
by the boundary conditions that the target imposes to the incident light field [||. Quantum- mechanically, Dicke [l[ 
has shown that, when two-level atoms are confined inside a volume much smaller than a radiation wavelength, the 
emission can be superradiant or subradiant. How the classical and quantum pictures are linked has been demonstrated 
at the example of a sample of weakly excited atoms @. More precisely, when a single atom out of N is excited, 
the Dicke symmetric state of maximum cooperation radiates superradiantly, i.e. at a decay rate proportional to N . 
Cooperative effects related to the superradiant and directional emission by an extended ensemble of atoms in a timed 
. . ■ symmetric Dicke state have been observed in the radiation pressure force acting on a large cloud of atoms driven by 
^ \ a resonant radiation field. Depending on the detuning of the incident radiation frequency from atomic resonance, the 
k>( • radiation pressure force may be either drastically reduced due to both increased forward scattering and a reduced 
scattering cross section [1, 0] , or even enhanced if the cooperative Mie scattering dominates over superradiance Q . 

In this paper we revisit the scattering by N atoms driven by a constant uniform radiation field and emitting 
radiation into free space. Our description assumes several approximations: 1) weak excitation of the atomic ensemble 
(one atom out of N is excited); 2) Markov or 'rapid transit' approximation Q (photon time of flight through the cloud 
much shorter than atomic decay time); 3) atoms frozen (zero temperature) and motionless; 4) neglected dipole-dipole 
interactions and collisions; 5) no resonant interactions, i.e. no resonance fluorescence, Van der Waals interactions etc. 
Approximation 3) excludes several cooperative effects related to atomic recoil motion (e.g. collective atom recoil lasing 
[9l [io| or matter wave superradiance [HI, [13), and thus neglects stimulated scattering processes along preferential 
directions, as for instance end-fire modes in Bose-Einstein condensates [l3[ or optical cavity modes |l4j . Neglecting 
atomic interactions (approximation 4) is justified assuming atomic distances much larger than an optical wavelength. 
The yet interesting opposite regime (atoms closer than an optical wavelength) requires the solution of the atomic 
equations with an exponential interaction kernel (see Eq.(j3j) and will be discussed in a future publication. 

We here determine the stationary solution for a spherical Gaussian distribution, as well as relevant macroscopic 
quantities such as medium polarization, scattered radiation intensity and radiation pressure force. Our solution 
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is based on a solution of the eigenvalue problem in the smooth density approximation. For large size samples, a 
continuous spectrum limit allows to obtain explicit analytical expressions for such quantities, expressed as a function 
of experimentally controllable parameters such as frequency and power of the driving field, optical thickness and size 
of the atomic sample. 

The paper is organized as follow: in Sec. II the scattering problem is expressed in terms of multi-particle coupled 
differential equations. In Sec. Ill the continuous density approximation is introduced and the solution for the atomic 
field is found in terms of the discrete eigenvalues of the interaction operator. Then, some macroscopic quantities 
of importance are calculated in Sec. IV. A continuous spectrum approximation for large atomic clouds performed in 
Sec.V appears to be particularly suitable to evaluate macroscopic quantities for arbitrarily large value of the atomic 
density. Sec. VI revises the radiation pressure force obtained assuming a symmetric timed Dicke atomic state [lj and, 
comparing this state to the exact solution, we find that phase shifts induced by the atomic cloud's refractive index are 
at the origin of important corrections for the expectation values of macroscopic quantities. Finally, numerical results 
and conclusions are presented in Sec. VII. 



II. COOPERATIVE SCATTERING PROBLEM 



A system of two- level (g and e) atoms with resonant frequency uj a and position tj, driven by an uniform laser 
beam with electric field amplitude Eq, frequency ujq and wave vector ko = (cJo/c)e 2 , is described by the interaction 
Hamiltonian: 
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where Qo = dEo/h is the Rabi frequency of the incident laser field, & j is the lowering operator for atom j, Ok is the 
photon annihilation operator and g^ = {cpLo^^heoVph) 1 ^ 2 is the single-photon Rabi frequency, where d is the electric- 
dipole transition matrix element and V p h is the photon volume. A special case is when a single photon is present in 
the mode k, as was extensively investigated in Refs.pl. [l5l fl6|. The total system (atoms+photons) is assumed to be 
in a state of the form |Q|: 
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where Ao = luq — Lu a - The first term corresponds to the initial ground state without photons, the sum in the second 
term are the states where a single atom has been excited by the classical field. The third term corresponds to the 
atom returned to the ground state having emitted a photon in the mode k, whereas the last one corresponds to the 
presence of two excited atoms and one virtual photon with 'negative' energy. It is due to the counter-rotating terms 
in the Hamiltonian (fTJ) and disappears when the Rotating Wave Approximation (RWA) is made. In the linear regime 
(i.e. a ~ 1) and in the Markov approximation (valid if the decay time is larger than the photon time-of-flight through 
the atomic cloud), the problem reduces to the following differential equation [l9j . 
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where T — Vphgy^/'KC is the single-atom spontaneous decay rate. The kernel in the last term of Eq.([3]) has a real 
component, — (T/2) YJ m ^ •[sin(/?j TO )//?j TO ] (where pj m — ^o| r j — r m|), describing the collective atomic decay, and an 
imaginary component, i{T /2) J2 m =ij[ C0S (Pj m ) I 'pjm]> describing the collective Lamb shift due to short range interaction 
between atoms, induced by the electromagnetic field [l9l - [2l| . The latter becomes significant when the number of atoms 
in a cubic optical wavelength, n\ 3 , is larger than unity, in which case the contribution from the virtual processes 
described by the counter-rotating terms in the Hamiltonian becomes relevant. Hence, for sufficiently dilute system, 
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such that N <C c 3 where a = fcoo\R and an is the cloud size, the collective phase shift arising from the imaginary part 
of the kernel in Eq. can be disregarded [22j and the scattering problem reduces to 

with initial condition /3j(0) = 0, for j = 1, . . . , N. Notice that Eq.Q deduced in the quantum mechanical description 
may be also obtained classically when the two- level atoms are treated as weakly excited classical harmonic oscillators 
0,0. For this reason the solution of Eq.Q (or of the approximated version ^) has a wider interest for the general 
problem of collective radiation scattering. 

As for the radiation field operator o,k, it evolves according the Heisenberg equation 

^ = ±[a k ,H] = -ig k e^-^ f] a m e~^, (5) 

m— 1 

where the fast oscillating term proportional to exp[i(a>fc + oo a )t] has been neglected. 



III. CONTINUOUS DENSITY APPROXIMATION 



In light scattering experiments, disorder (or granularity) plays a role when the number of atoms projected onto a 
cross section perpendicular to the incident beam is small enough so that a light mode focused down to the diffraction 
limit (that is ~ A 2 ) would be able to resolve and count the atoms. In other words, the stochastic fluctuations induced 
by the random positions of the atoms can be neglected when the total number of atoms N is larger than the number 
of modes ~ a 2 that fit into the cloud's cross section, i.e. when the optical density is bo = 3N/a 2 ^> 1. Under this 
hypothesis, the particles can be described by a smooth density n(r), and their probability to be excited by a field 
/3(r, t), so that Eq.Q turns into 

i/)(r, t ) = (iAo - «r,f) - 4<*" -\j d ,„ (r0 5i^L-il /V ,^ (6) 

In what follows we will consider only spherically symmetric distributions n(r). Because of the linearity of Eq.©, it is 
convenient to introduce an eigenbasis of the coupling operator. The functions j n (r)Y nm (Q, cj>), with j n s the spherical 
Bessel functions and Y nm {9, <p)s the spherical harmonics, appear as a natural choice considering the following identity 
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MH(/ ' n|1 ' n = 47r E E Ukor)Y^ n (0,4)Y nm (e , ,^)j n {k o ^). (7) 



fcn|r — r . 

1 1 n— o m=—n 

In particular, the choice of the spherical harmonics guarantees the orthogonality of the basis, since 
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Therefore, assuming the following decomposition for the field, 

oo n 

/8(r,t)=£ <Xnm(t)jn{k O r)Y nm (0,(t>), (9) 

n— m=—n 

the projection of Eq.© along the eigenmodes leads to 
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jn(k r) = --O nm , (10) 



where A„ are the eigenvalues associated to modes n 

A„=47r/ r 2 n(r)j 2 (k r)dr (11) 
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whereas £l nm corresponds to the projection of the incident wave on mode (n, m) 
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n nm = n / d0 / d9 sin 9 ¥^(9, 0)e lk ° rcos 6 = 2n 8 m oV<^ + l)i n j n (k r). (12) 



Assuming the cloud is initially unexcited, i.e. a mn (0) — 0, only spherically symmetric components with m = will 
grow so that, defining a n (t) = a n o(t), Eg. (ITU)) reduces to 

r 
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a n = -i n+1 y/n(2n+l)Slo (13) 
Eg.([T^l) straightforwardly integrates and, inserted in Eq.©, leads to the following expression for the excitation field 
M r , 9,t) = ^ Y * n ( 2n + 1 )jn(kor)Pn(cos6) r _ (r/2)(1+A „ )t j ^ {u) 

n— v ' 

where the scaled detuning 8 = Aq/T was introduced. Hence, each mode n relaxes toward the steady-state with a 
characteristic time r n = l/r(l + A n ): the first modes relax very quickly since A„ is proportional to N, yet for the 
highest modes, t„ ~ T -1 , even if their macroscopic contribution is usually small. Eventually, for times much longer 
than the single-atom decay time the field tends toward a stationary state fully characterized by the spectrum 

*<'•« - Tt»iT+l) M, * r)p " lm *>- (15) 

n— 

Notice that the set of eigenvalues ffTT]) is complete since from the identity J2 n>0 {2n + l)j^(z) = 1, it follows that 

V(2n + 1)A„ = 4tt / r 2 n(r)dr = N, (16) 
n=0 Jo 

which corresponds to the trace of the coupling operator. 

IV. MACROSCOPIC QUANTITIES 

The description of the field /3 s (r, 9) in terms of spectrum also provides expressions for any macroscopic quantities, 
the most relevant of which are calculated here below. These formulae will be specialized to Gaussian clouds in the 
subsequent section. 

A. Average amplitude and probability of excitation 

The average 'phased' probability of the timed Dicke state [lj| and the excitation probability are respectively 

(/3 s e- lk °- r ) = — / d<9sin<2 / drr 2 n(r)p s (r, 0) e - tk " rcose (17) 
N Jo Jo 

(|/? s | 2 ) = — / dflsinf? / drr 2 n(r)\P s (r,9)\ 2 . (18) 
Inserting Eq. p3|) and using the identities 
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we obtain 
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B. Scattered field 



The electric field radiated by the excited atoms reads, in the smooth density limit (see Appendix A), 
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Using the stationary solution (1151) for a spherically symmetric distribution and using the integral 

/ d9' sm9' P n (cos9')J (k r' sm9 sm9')e- lk ° r ' cosf> cosf> ' = 2r n j n {k r')P n {cos9) 
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we obtain 



(22) 



(23) 



(24) 



where we used the relation V — d 2 kg / (27r& ) . In the forward direction (9 = 0) the scattered field is proportional 
to the 'phased' probability amplitude (|20|) . The expression ([24]) provides the angular distribution of the scattered 
radiation field. We can also calculate the scattered intensity as 
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Passing to the continuous limit and using Eqs. ([2~Tj) . (f2"2")l and (IM1) . we obtain: 
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where Jo — ceoE 2 . The scattered intensity is the sum of the incoherent contribution, proportional to N (since A„ oc N, 
see Eq.([TI|) and isotropic, and the superradiant contribution, proportional to N 2 and directed (for extended clouds) 
mainly in the forward direction. Integrating over the solid angle, the total scattered power is 



P s = 2irr 2 [ d9 sin 6I s {r, 9) 
Jo 

where we used the second of the identities (|T9")) . 
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C. Radiation pressure force 

The radiation force operator acting on the jth atom is calculated from Eq.flT]) as Fj = —V r .H = F a j + F e j where 

1 

F aj = ^k ^{a J e i ( At ' t - kor ^-h.c.} (28) 
F ej = ^Ek5fc{4^e l(Wfc "" a)t " lk ' rj -<7]a k e- l(Wfc -^ )t+4k -^} (29) 

k 

are the forces due respectively to the absorption and emission processes. In Eq. ([29|) we have neglected the 
counter-rotating terms proportional to exp[±i(wfc + uj a )t\. We are here interested in the average radiation force 
F = (1/iV) J2j Fj, which stands for the force acting on the center-of-mass of the atomic cloud along the direction 
of the incident field, ko = k$e z . This average force is relatively easy to measure by time-of flight techniques in 
cold atomic clouds released, for instance, from magneto-optical traps (MOTs) and has recently revealed cooperative 
effects in the scattering by extended atomic samples It may provide a convenient measurement (aside from the 
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scattered radiation) of the effects that cooperative scattering imprints on the atoms. The average absorption force 
along the z-axis, resulting from the recoil received upon absorption of a photon from the incident laser, is 

N 

F a = = ^hk n [<T,e jAot - lk ^ - h.c] . (30) 

The second contribution F e = (1/N) F e j results from the emission of a photon into any direction k. Inserting 
from Eq.([5]) into Eq. (|29p and approximating the sum over the modes k by an integral, it is possible to obtain, in a 
way similar as done for the radiation field operator Es of Eq. (|A5|) . the following expression for the average emission 
force along the z-axis 0: 

F e = - j^ j d<t> f dflsinflcostf ^ [e-^^-'^a^aj + h.c] . (31) 

j^rn—l 

Evaluating their expectation values on the state @ (neglecting virtual photon contributions), the emission force in 
the discrete model is 

(Fe) = p dcf>r d9 S m6 cos6 [ft/^e-^-^+c.c 

= *W E j^7^ii(fco|r,-r ro |)(«-c.c). (32) 

where we used the identity (|B1[) and ji(z) is the first order spherical Bessel function. Then, passing to the continuous 
distribution limit, Eqs. (|3"0")) and (|3"Tj) are approximated by 

(F a ) = -hk n lm(l3e~ lk °- r ) (33) 

(F e ) = -Mo-^ J d(f> J d0sin0cos0 J dvn(r) J dr'n(r') [/3(r)/3*(r')e^ k (l "- r,) + c.c] . (34) 

The absorption stationary force is readily obtained from Eq.(|20|). 
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whereas a longer calculation, reported in Appendix B, yields the emission stationary force, 

/p\ *i> 2n ° V- (" + VKK+iW? + (1 + A n )(l + A n+1 )] 

{Fe) = °rWh [4<5 2 + (l + A„P][4^ + (l + A„ +1 ) 2 ] ' (36) 

As expected, the absorption force pushes the atomic cloud in the direction of the driving field, whereas the emission 
force is oppositely directed and is proportional to N. Both forces depend on N also through the collective decay rate 
T(l + A„). We observe that the absorption force ((33)) is linked to the scattered power (|2"Tj) by the relation 

N(F a ) = Att (^-) , (37) 



i.e. the absorption force is proportional to the scattered power per atom, Pg/N. 

Fig Q] compares the (a) absorption, (b) emission and (c) total forces vs. N for S — 10 and a Gaussian density profile 
with cr = 8 (where a = kocru), calculated using the series (|35|) and (|36|) (circles) and iV-body simulations (see EqHJ 
(dots). The eigenvalues A„ for the Gaussian density profile are given by Eq. (|3"5|) . The forces are reported as a ratio 
between the cooperative force and the single-atom force, F\ — hk^TVl^/ (4A 2 , + T 2 ). The error bars correspond to the 
standard deviation of the observables over 8 realizations. FigQJl shows the total force vs. N for a different choice 
of parameters, a — 5 and 6 — 200, for which the force exhibits a maximum as a function of N [8J. We observe a 
good agreement between the analytical solution (|3"5j) and the TV-body simulations at large N, when effects due to the 
discreteness of the system become negligible. 
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FIG. 1: (a) Absorption, (b) emission and (c) total forces vs N for 5 = 10 and a Gaussian cloud with a = 8. (d) Total force 
for <t = 5 and <5 = 200. Forces are relative to the single-atom force F\. The plain curves refer to the series (|35[1 . (|36|l and 
their sum, the circles to the analytical expressions (|46[) , (|36|) and their sum, the dots to iV-body simulations (see Eq(4]| and the 
dash-dotted lines to the STD state. The error bars correspond to the standard deviation of the observables over 8 realizations. 



We recall that our expressions have been obtained in the continuous density approximation, i.e. assuming a sample 
with large optical thickness, i.e. N ^> a 2 , but sufficiently dilute to neglect the collective Lamb shift, i.e. with a small 
number of atoms in a cubic optical wavelength volume, N <C <J 3 . These two conditions imply that the cloud size should 
be much larger than the optical wavelength, a 1. Nevertheless, our results remain valid also in the limit n(r) — > 
where, in particular, the standard radiation pressure force, F\, is recovered. However, at low densities, fluctuations 
will strongly affect individual measurements of the radiation pressure force and in this regime our expressions represent 
only the average expectation values. 



V. LARGE GAUSSIAN CLOUDS 



Large clouds (er 3> 1) behaves fundamentally differently from small clouds (a < 1), as can be deduced from their 
spectrum. Let us for example consider the case of Gaussian clouds, with density n(r) — (N/(2n) 3 / 2 a^ i ) exp(— r 2 /2cr|.). 
The spectrum then reads 

A " = N \^r P 2e ' p2/2 ^P)^P = N ] J^e-° 2 I n+1/2 (a 2 l (38) 

where a = fcocr^ is the scaled size of the cloud, and I n (x) the n-th modified Bessel function. Recently the spectrum A n 
for an uniform spherical cloud has been calculated by A. A. Svidzinsky et al. [1, also for the exponential kernel of 
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Eq.Q. However, a Gaussian distribution is certainly more realistic for experiments with cold dense atomic ensembles. 
Generally, the spectrum of small clouds (cr < 1) is composed of a few significant eigenmodes, whereas for a large, all 
the eigenmodes for n < a are significant and the spectrum can be treated as a continuum. In particular, in the latter 
case the A^ can be approximated, for n < cr, by A„ ~ (N/2a 2 ) exp[— (n + l/2) 2 /2cr 2 ] (see e.g. Switching to a 

continuous treatment of the spectrum, we define rj = n + 1/2 and get 
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Remark that using these definitions, the completeness condition (|16l) is still preserved. The continuous spectrum limit 
allows for the evaluation of the sums in Eqs. (l2T)l) and (|2"TT) as continuous integrals 
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where we have set x = (bo/6) exp(— i] 2 /2a 2 ), with bo = 3iV/cr 2 the optical thickness. The above expressions integrate 
as 
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These formulae highlight the prominent role of the parameters bo and 5 in the high-density limit. In a similar way, 
we calculate the total scattered power as: 
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As expected, for small optical thickness (bo <C 1) the scattered power is incoherent, Ps ~ [/o^V/(47r/c 2 )]/(4(5 2 + 1). 
However, for large optical thickness it shows a logarithmic dependence on N . The superradiant character of the 
radiation is visible only observing the scattered intensity in the forward direction (see the second term of Eq. ([26|) ). 
but not in the total scattered power. 

The absorption force is deduced from Eqs. (|3"3"|) and (|4"3")l as 
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The emission force can be written, in the continuous spectrum approximation in the integral form 
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,(47) 



where we have set y = (77 + l)/cr. As can be observed in FigfTJ the continuous-spectrum approximation gives excellent 
results compared to the full series ([35]) and (|36p . In the limit a — > 00 and finite &o, Eq. (|4T|) would lead to 



(F e )^~(F a )+nk T(\(3 s \ 2 ), 



(48) 



which has a transparent interpretation: for a very large cloud, the atoms scatter radiation to forward direction and 
the recoil received by the atoms upon emission cancels out with the recoil received by absorbing a photon from the 
driving field. The net force remaining after the substraction is the non collective contribution to the emission. The net 
force is equal to the photon momentum hko times the emission rate r(|/3 s | 2 ). This emission rate depends indirectly 
on N and a through the enhanced superradiant decay, T(iV/4a 2 ), which decreases the emission rate when the optical 
thickness increases. A more accurate expression of the emission force valid for large but finite cloud size would require 
the exact evaluation of the integral in Eq. (|4T|) . 
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FIG. 2: Average excitation probability (|/? s | 2 } (left) and variance up — (|/3 S | 2 ) - |<^e~ lfc °- r >| 2 (right) vs AT. The plain curves 
refer to the analytical expressions (I43[) and (1441) . the circles to series expressions (I20[) and (1211) . the dots to ./V-body simulations 
(see Eq[4]) and the dash-dotted lines to the STD state. The error bars correspond to the standard deviation of the observables 
over 8 realizations. Note that for the STD state, ap = 0; for too small N, the approximations (|43[1 and (|44l) provide inconsistent 
results, that is negative a p. Simulations realized for a — 10 and 5 = 10. 



VI. SYMMETRIC TIMED DICKE STATE 



A particular ansatz used by Scully and coworkers @, H [l6| is the symmetric timed Dicke (STD) state, given by 

/3(r,i)=/3 T D(t)e lk °- r . (49) 
After integration over space of Eq.f6]), one obtains the following evolution equation 



dt 



iA - ~(l + N Soo ) 



(50) 



where Soo is the integrated structure factor of the cloud defined as 



, x = _L / n ( r ) dr / dr , n(r / ) aiii L V 7J 1 J t [Z e -»ko ■(■—-') = ±_ j d(j) j d9 sm9\(e l{k °- k > r }\ 2 . (51) 



„/^^N Sin ( fc 0| r - r, |) c -»k o .(r-r') = 1 

ko\r — r'l 47r 



2'K 



This ansatz is of particular interest since it allows to evidence e.g. the superradiant nature of the decay of such state 
when the pump is turned off . As for its steady-state, it reads [3, 0] 



Ptd — 



fi 



1 



r 2S + i(l + N Soo )' 
and for a large cloud with Gaussian distribution, since Soc « l/4cr 2 , we get 

1 

PTD r 28 + t(l + 6o/12)' 
Thus, the STD solution ([53)) approximates the exact result (j44|) only for &o "C 12, 

1 



TD 



r 



45 2 + (1 + bo/6) ' 



(52) 



(53) 



(54) 



Fig J2] shows the average excitation probability < | /3 S | 2 } (left) and its variance, er 2 , = (\(3 S \ 2 ) — \(f3 s e lfc o- r )| 2 vs. N 
for a — 10 and S = 10. The plain curves refer to the analytical expressions (|43j) and f|44[) . the circles to series 
expressions (|20p and (p?Tj) . the dots to A^-body simulations (see EqU]) and the dash-dotted lines to the STD state (|49l) . 
The error bars correspond to the standard deviation of the observables over different realizations. Note that for the 
STD state, op = 0; for too small N, the approximations (|43l) and (|44l) provide inconsistent results, that is negative 
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erg. We observe excellent agreement between the series and the analytical solutions, and a consistent reduction of 
the excitation probability decrease vs N with respect to the TDS prediction (dashed line in figHJ left). Also the 
fluctuations obtained from the A^-body simulations converge for large values of N toward the results obtained in the 
continuous density approximation (figlH left), showing the presence of a shot-noise contribution for small N. 
The radiation force for the STD ansatz @J is dHB 

(F) = (F a ) + (F e ) = hk [-SIqMPtd) ~ r|/3 T D| 2 A7oo] (55) 

where 

/• 2 77 



/oo = i /"(ty r d8 sin6 cos 9\{e^ k °-V r )\ 2 . 
4tt Jo Jo 



(56) 



Since for a spherically Gaussian distribution N(soo — /oo) m Nf^a 4 ) — bo/(24a 2 ), from Eq. (|52| the stationary 
radiation force is 0, 0, Hi 



n \ 2 l + N( SoD -f a 



" ^HTj 4^ + (l + Vl2) 2 - (58) 

As can be observed in FigfTJ the STD state yields a good agreement with the full-spectrum approach and the iV-body 
simulations only for small values of the optical thickness. 

Finally, for that state (l49l) the scattered radiation electric field, Eq. (|A5p , and intensity, (|25|) . become 

e s (r,6^,t) = ^^ e ^o(r-ct)^ TD(t)(e i(ko-k). r) (5Q) 



Aireor 



and 



( dkl 



2 



Is(r,M) = ceo ( ^ ) \PTD\ 2 [N + N'\(e^-^)\ 2 \- (60) 

In particular, for a Gaussian distribution, (exp[i(ko — k) • r]) = exp[-cr 2 (l - cos 0)/2]. 

Fig [3] shows the phase of the excitation amplitude f3 s (left) and the excitation probability, \(3 S \ 2 , (center) in the 
(x, z) plane (y = 0), calculated from the exact solution (fT5|) . for = 10 4 , cr = 10 and 5 = 10. FigJ3] (right) shows the 
contribution to the electric field radiation in the same plane, weighted by the local atomic density. The simplification 
of the STD, as compared to the exact solution (|15p . resides in the assumption that all atoms are equally excited 
and oscillate in phase. According to the exact calculation, the atomic dipoles appear to be in phase only in the core 
of the cloud (see Fig[31 left), but this phase profile has strong distortion away from it. This phenomenon is all the 
more important as the atoms are much more excited in the peripheral region than in the core (see Fig|31 center). In 
particular, even when this excitation probability is weighted by the particle density, two areas at the cloud entrance 
and exit contribute significantly to the radiation electric field (IM1) (see peaked structures in FigEl right). For a timed 
Dicke state, both the phase profile and average excitation remain fiat throughout the cloud. From a macroscopic 
point of view, the STD state neglects phase shifts imprinted into the pump beam by the cloud's reflective index. This 
can be seen in Fig. Q|a), which compares the phase of the STD state, Ptd exp(ik • r), (linear curve, no phase shift) 
and /3 s (r) (additional phase-shift) along the optical axis across the cloud. Fig. 0Jb) shows the pump beam phase shift 
after transmission through the atomic cloud as a function of atom number. This phase shift is at the origin of the 
deviation between the radiation pressure forces calculated for the STD state and the exact solution. The pump beam 
phase shift leads to a reduction of the absorption and the emission forces. This can be understood as destructive 
interference of forward radiation emitted from different atoms, located at the same plane z — zq but different x or y. 

As for the emitted wave, it is concentrated in the forward direction (see Fig|5j left), and there is no backscattering 
(not shown here). The wavefront phase does not exhibit significant distortion in the central region of the radiated 
beam (see FigEl right). 



VII. CONCLUSIONS 



In summary, we studied collective scattering from a dense and large atomic cloud with gaussian density profile in 
terms of the eigenvalues of the interaction operator. This enabled us to calculate the state generated by the interaction 
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FIG. 3: Phase of the atomic radiation ang(/3 s ) (left), Excitation probability |/3 S | 2 of the atoms in the (x,z) plane (y = 0) 
(center) and contribution to the radiation of the electric field in the same plane, which corresponds to the level of excitation of 
the atoms weighted by their local density (right). Simulations realized using the analytic expressions (|15[) , for TV = 10 4 , a = 10 
and 8 = 10. 




FIG. 4: Left: Phase shift along the optical axis for the exact solution (|15[) (plain line) and a STD state ((49} (dash-dotted 
line). Right: Phase shift of the pump beam after transmission through the cloud as a function of TV evaluated along the axis , 
x = y = 0, at k z = 20a, for TV = 10 4 , a = 10 and 8 = 10. 

with a laser beam. We found that this state considerably deviates from the 'timed Dicke' state at high optical densities. 
In order to characterize this state, we calculated the phase front of the excited atomic dipole moments and the force 
due to cooperative scattering by the atomic cloud. 

Our approach consisted in expanding the solution on spherical harmonics, which form a complete orthogonal basis 
of the Hcrmitian interaction matrix. Under the assumption that the system is sufficiently dilute to neglect short-range 
effects due to dipole-dipole interactions, the continuous spectrum approximation allows to transform infinite series 
into solvable integrals, and eventually to derive analytical expressions for the most relevant observables, such as the 




FIG. 5: Transverse intensity profile (left) of the wave emitted by the cloud, at z — 3<r and its wavefront (right). Simulations 
realized using the analytic expression (|24p . for TV = 10 4 , a = 10 and 8 — 10. 
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scattered intensity and radiation pressure force. These analytical expressions show good agreement with the numerical 
solution of the A-body problem, and highlight the dependence on the accessible experimental parameters such as 
optical thickness, atomic cloud size and laser frequency. 

The analytical solution appears particularly useful for studying the thermodynamic limit when N —> oo and 
V — > oo with N/V fixed, until collisions or non resonant interactions come into play. The thermodynamic limit is 
hardly accessible to TV-body simulations, since the latter are highly CPU consuming. 

In contrast, the eigenvalue approach opens the possibility to study the fascinating link between microscopic and 
macroscopic domains of light scattering, and in particular between single-atom scattering and Mie scattering for 
extended continuous samples @. 

For large optical thickness, the refraction index of the cloud acts back on the driving field and shifts its phase. 
For this reason, the solution for large optical thickness shows appreciable deviations from that obtained assuming a 
symmetric timed Dicke state for the atomic sample, since the latter corresponds to a completely degenerate eigenvalue 
spectrum. The exact solution (| 15[) takes into account the induced phase shift and reduces to the symmetric timed 
Dicke state in the limit of relatively small optical thickness. 

An important further development of the present study should be to understand how the dipole-dipole interactions 
contribute to the observed cooperative effects, completing in this way the cooperative scattering description for highly 
compressed and dense atomic clouds. 
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Appendix A: Derivation of the scattered electric field 



The scattered radiation is provided by the positive-frequency part of the electric field operator 

E s (r, t)=J2 £kMty k ' r ~ iUkt (Al) 

k 

where £ k = (fku k /2eQVph) 1 ^ 2 ■ Integrating Eq.© with a^(0) = 0, inserting it in Eq. (|Al[) and approximating the sum 
over the modes k by an integral, we obtain 

E s (r,t) = -i^g J2 J o dt'a m (t - t'y^' J dk£ k g k e tk < r - r ^- tckt ' . (A2) 

Introducing spherical coordinates, dk = dkk 2 d(f>d9 sin 9, and integrating the angular part, Eq. (|A2p becomes 

E s (r,t) = V f dt'a m (t - t'Y^' ^ dkk£ k g k { e -**(f-k-r«|/c) _ e -**(f +|r-r m |/c)| (A3) 

m =l I 1 " ~~ Vm ' JO JO 1 ' 

Assuming the radiation spectrum centered around k w ko, we approximate k£ k g k ~ ko£ ko g ko . Then, extending the 
lower limit of integration of k to — oo, we obtain for t < |r — r m |/c [!, [23[ 

Mr,*)«- ^EyT7j^-|r-r m |/c). (A4) 

m— 1 1 1 

Neglecting the radiation retard in the limit t ^> cr^/c, where an is the cloud size, and approximating in the far field 
limit |r — r m | « r — n • r m where n = r/r, Eq. (IA4l) becomes 

dk 2 N 
E s {r,6,4>,t) « _i!L e ft»M)4^t V a m (t)e-*—, (A5) 
47reor L — ' 



m— 1 



where k s = fco (sin 9 cos cj>, sin 9 sin <f>, cos ff) . When applied on the state ©, neglecting virtual transitions, it yields 
Es\ty) — £s\gi ■ ■ -9n} an d where in the continuous density approximation will be given by Eq. 
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Appendix B: Derivation of Eq.(|36[) 



The angular integral in Eq. (j34|) is 

/•27T pTC 



f C 7 — z' 

/ d<t> / desm6cos6e- %k ( - r - r) =4ni- jrjMr - r'|). (Bl) 

Jo Jo |r - r 



Since 



d (z- z') 

-^-jo(k \r - r'|) = -/son 7T-?i( fc o|r - r'|), 

02; |r — r| 



where jo(x) = sm(x)/x, Eq. (|34f can be written as 

(F e ) = -i^ | drn(r) {/?«^y / dr'n(r')io(fco|r - r'|)/3*(r') - c.c.j . (B2) 
Using the expansions ([7]) and © and Eqs.© and (fTTj) , we obtain 

„ 00 n 

/ dr'n(r')io(fc |r - r'|)/T(r') = J2J2 <m*njn(k r)Y* m (8 , 0). (B3) 

w n , ~ 



In spherical coordinates 



n—0 m=—n 



d „ 3 sin 9 d 

— = cos 9— — h 



dz Or r dcosO 

and 

,2 



(# e ) = -i^/dmCr) £ < m A„(co S ^ + ^^)i„(M^ m (M)-c.c. (B4) 

^ n—0 m— — n. \ \ / ✓ ^ 

Still using Eq.©, 

{Pe) = ~^r£ EEE A » drr 2 n(r) Jp (fc r) / d0 / dflsinflx 

x | d(fc r) cosgl Pg^0) y nm( g .0) H ^-sm gr w (g,^) d( C os(9) h(<W*nm-c.c.). (B5) 

Assuming as before 

and since 



/ 2n + 1 

r n0 (M) = ^/^^p„(cos( 



Eq. (IB5|) becomes 

(Fe) = 
X 

Since 



Vy'A n V'(2n + l)(2p+l) / drr 2 n(r)j p (k r) d9 sin 6 x 
P =0 n=0 ^° ^° 

«9j„(fc r) cos ^ Fp(cos ^p ji(cos0) + S in 2 6 P p (cos g) ^ (C °" f } } i (o p «; - c.c.) . (B6) 

f^O^ CI COS (/ ) J 



hkoT 

' AN 

p=0 n=0 



d(k r) 



J d9 sin 6* cos 9 F p (cos 6*) P„(cos 9) = J dxxP n (x)P p (x) 
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and 



using the identities 



d9 shr 3 6 PJ cost 



I dP n (cos 6 
<9(cos#) 



cLc(l-a; 2 )Pp(a;) — P„(s), 
i da; 



and 



we obtain 



and 



(2p+l)xP p {x) = (p+l)P p+1 (x) + pP p ^(x), 

dP^x)_ = n 

ax 



dxP n (x)P p (x) = - — — S np , 
l 2n + 1 



dxxP„ (x)P n (x) — — r 

_ x * U pU (2p+l)(2n+l) 



f 1 d 

J dx(x 2 -l)P p (x)—P n (x) 



2)i 



(2n+l)(2p + l) 
By substituting these expressions in Eq. (|B6[) we obtain 



{(j> + l)*n,p+l + P S p,n+l} 

{p$ P ,n+i - {n + l)S n , p+ i} . 



2iV 
x Up+1) 



P =o„=o v/(2n+l)(2p+l) Jo 

9j P +i(k r) O jp+i( k or) 



drr 2 n(r)jp(kor) x 



+ (p + 2) 



Since 



and 



5(/c r) ' k r 

djn+l(z) 



S n ,p+i + P 



dj n (k r) j n (k r) 
n 



d(k r) k r 



+ 1 



dz 

djn(z) 



■ i \ 71 + 2 ■ n 

Jn{Z) ]n+l{z) 



using the definition (JTTJ, 



^^i(a p a; -c.c.) 



p=0 n=0 

Eliminating one of the two sums, we obtain: 



A /(2n + l)(2p+l) 



{(p + l)5„,p+i - _p(5 p ,„ + i} . 



E(n + l)A„A„+i 
— ^======i(a„a„ +1 - c.c. i 



4vr7V ^ ^(2n+l)(2n + 3) 



In the stationary case, 



and the stationary force is 



do i n y/4n(2n+l) 

T 2<y + i(i + a„) 



/p\ = -ft*. 2n " V (" + l)AnA»+i[4^ 2 + (1 + A»)(l + A n+ i)] 
{ e) ° rA ™o [4<5 2 + (l + A„) 2 ][45 2 + (l + A„ +1 ) 2 ] ' 



(B7) 



(B8) 



(B9) 
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